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ABSTRACT 



Long-range scalar forces with gravitational strength can be ruled out by accurate tests 
of scalar-tensor theories in the solar system. However, such tests are based on the motion of 
celestial bodies made of standard baryonic particles, which leaves room for possible scalar 
interactions in the dark sector of the universe if a new scalar degree of freedom couples 
only to dark matter particles. In particular, an interaction between Cold Dark Matter (CDM) 
and a classical scalar field playing the role of the cosmic dark energy (DE) might provide 
such long-range dark interactions without conflicting with solar system bounds. Although 
presently available observations allow to constrain such interactions to a few percent of the 
gravitational strength, some recent studies have shown that if CDM is composed by two 
different particle species having opposite couplings to the DE field, such tight constraints 
can be considerably relaxed, allowing for long-range scalar forces of order gravity without 
significantly affecting observations both at the background and at the linear perturbations 
level. In the present work, we extend the investigation of such Multiple Dark Matter scenarios 
to the nonlinear regime of structure formation, by presenting the first N-body simulations 
ever performed for these cosmologies. Our results highlight some characteristic footprints 
of long-range scalar forces that arise only in the nonlinear regime for specific models that 
would be otherwise practically indistinguishable from the standard ACDM scenario both 
in the background and in the growth of linear density perturbations. Among these effects, 
the formation of "mirror" cosmic structures in the two CDM species, the suppression of 
the nonlinear matter power spectrum at fc > 1 h/Mpc, and the fragmentation of collapsed 
halos, represent peculiar features that might provide a direct way to constrain this class 
cosmological models. 
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1 INTRODUCTION 

While the existence of two distinct forms of gravitating energy in 
the Universe which do not interact with the elecromagnetic field 
- named Dark Matter and Dark Energy - is now supported by a 
large number of independent and complementary observations (see 
e.g. |Bertone et aT1|2005l |Bartelmann|[20l0l |Astier & Pain1|20~i2l 
and references therein), their fundamental nature is still completely 
unknown and represents the main focus of most of the present 
research efforts in the fields of theoretical and observational cos- 
mology and astroparticle physics. Most of our knowledge about 
the phenomenology of Dark Energy (DE) and Cold Dark Matter 
(CDM) comes from the direct observation of their gravitational ef- 
fects on the evolution of cosmic structures and on the global evolu- 
tion of the Universe as a whole (see e.g. Riess et al. 1998 , Schmidt 
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etal.|19 98 , Perl mutter et al. 1999 , Perc ival et al.|2 001 , Astie r et al.| 
20061 rGiannanto nio et al.]20 08, Guzzo et al. 2008, Fu etal.|2008[ 



Mantz et al . 2010, Schrab back et al.|[20Tol [Percival et al.||20101 
Sherwin et al. 2011). However, such effects do not allow to drive 
clear conclusions about the fundamental constituents of these two 
dark components, and in particular do not allow to resolve a possi- 
ble internal complexity of each of these fluids. Furthermore, all the 
experimental and observational initiatives aimed at a direct detec- 
tion of CDM particles in the laboratory have so far failed to provide 
an unambiguous and statistically significant identification of any of 
the proposed CDM particle candidates (for an excellent recent re- 
view see e.g. Bergstrom 20121. 

Such a deep uncertainty on the fundamental nature of these 
two cosmic fluids leaves room for theoretical speculations about a 
possible hidden complexity of the dark sector, as e.g. a dynamical 
nature of DE based on the evolution of an ubiquitous scalar field 
as in the case of Quintessence (Wetterich 1988, Ratra & Peebles 
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[7988] > or k-essence (fArmendariz- Picon et al.|2001| l, a direct inter- 
action between DE and CDM (see e.g. Wetterich 1995, Amendola 



2000, Farrar & Peebles 2004 , Pettorino & Baccigalupi 2008 , Baldi 
2011b|l, a mult i-particle composition of the CDM fluid (see e.g. 



Khlop ovfT995l |Brookfield etaLp OOS , Baldi |2012b||Maccio' et alT] 



2012), or the presence of long-range scalar or vector interactions 
between CDM particles (see e.g. Hellwing et al.|2010{ Keselman] 
[eFaT1[20T0l |Loeb & Wemer]|2011[ |Aarssen et al.||2012| >. The irT 
vestigation of such theoretical scenarios in the last years allowed 
to put tight constraints on their characteristic parameters and to 
strongly limit the range of viable theoretical extensions of the stan- 
dard ACDM cosmological model. In particular, models featuring 
a direct interaction between a DE scalar field and a single fam- 
ily of CDM particles have been widely studied for what concerns 
their impact on the background expansion history of the universe, 
as well as for their effects on the growth of density perturbations in 
the linear (e.g. |Bean et al.|2008| |La Vacca et al.|20091|Xia|2009> 



and nonlinear ( |Macci6 et al.|2004[|Baldi et al.|2010[|Baldi|2011a| 
|Li&Barrow|2011| > regimes of structure formation, which allowed 
to put relatively tight constraints on the coupling constant between 
DE and CDM and consequently on the strength of its associated 
long-range scalar force. 



Nevertheless, the recent works of Brookfield et al. (2008J and 
|Baldi|p012b} have shown how in the presence of a higher level 
of internal complexity of the CDM sector, as e.g. in the case of 
multiple CDM particle species with different couplings to a DE 
scalar field, the background and linear constraints on such dark 
scalar forces might be significantly relaxed, even for scenarios that 
do not introduce any additional parameters with respect to stan- 
dard models of interacting DE. In the present work, we will extend 
the investigation of these coupled DE models with Multiple Dark 
Matter families to the nonlinear regime of structure formation by 
means of a series of intermediate-resolution N-body simulations, 
with the aim to identify possible characteristic footprints capable 
to constrain the parameter space of such models beyond the very 
loose bounds provided by the background expansion history and 
by the linear growth of structures. As our study will show, although 
nonlinear structures result substantially affected in these models for 
sufficiently large values of the DE-CDM coupling constant, a long- 
range scalar interaction with a strength comparable to that of stan- 
dard gravity cannot be easily ruled out even in the nonlinear regime. 
Such results open new possible directions for phenomenological 
extensions of the standard cosmological model, allowing for new 
long-range interactions of gravitational strength. The present work 
therefore represents a rather qualitative investigation of the effects 
of multiple interactions between DE and CDM in the nonlinear 
regime, and is mainly aimed at a broad sampling of the parame- 
ter space of these scenarios which could provide a useful guideline 
to drive future numerical investigations towards selected parameter 
combinations. 



Our paper is organized as follows. In Section|2]we will intro- 
duce the specific cosmological models with Multiple Dark Matter 
families interacting with DE. In Section[3]we will describe our set 
of numerical simulations and in SectionRfwe will present the main 
results of our analysis, concerning the shape of large-scale struc- 
tures and the nonlinear matter power spectrum in such cosmologi- 
cal models. Finally, in Section|5]we will draw our conclusions. 



2 MULTIPLE DARK MATTER MODELS 

If the observed accelerated expansion of the Universe is sourced by 
some dynamical scalar field such as Quintessence (|Wetterich|1988) 
|Ratra&Peebles|1988| > or k-essence (Armendariz-Picon et al. 2001 ) 
rather than by a simple cosmological constant, it is natural to con- 
sider the possible direct interactions that such DE component might 
have with the other fields of the dark and of the visible sector. Any 
such interaction would then mediate a new scalar force between 
massive particles, and in most viable scenarios this force is long- 
ranged since the DE field has to be sufficiently light in order to ex- 
plain the energy scale of DE, and its propagation length is therefore 
correspondingly large, generically comparable to the cosmic hori- 
zon (for a general review of interacting DE models, see Chapter 1 
of |Amendola et al.|2012) . 

Universal interactions between a cosmic scalar and all mas- 
sive particle fields in the universe constitute the broad family of 
scalar-tensor theories of gravity and Extended Quintessence mod- 
els <|Brans & Dicke|[T96T] [ Appleby & Welrer||2010| |Baccigalupi| 
|et al.|2000||Perrotta et al.|2000||Pettorino et al.|2005| >, and their in- 
teraction strength is very tightly constrained by solar system tests 
of General Relativity (see e.g. |Bertotti et al.|2003[|Will|2005> un- 
less some screening mechanism - as e.g. the Cham eleon {Khoury 



& Weltman 2004), the Vainshtein ( Vainshtein 1972 



2002), or the Symmetron ( Hinterbichler & Khoury 



[Deffayet et al. 



2010 1 effects - 



hides the scalar force in our local environment. However, such tight 
constraints do not apply to a selective interaction between a DE 
scalar field and the dark matter sector, as first proposed by Damour 



DE scenarios (see e.g. Wetterich 


1995; 


Amendola 2000 2004f|Pet-| 


torino & Baccigalupi||2008 Ba 


di||20 


Llb| 2012aJ that have been 



very widely investigated in the last decade as they could provide 
dynamical solutions to the fine-tuning problems of the cosmolog- 
ical constant. In these models, DE and CDM (which is assumed 
to be made of a single species of weakly interacting and massive 
fundamental particles) interact through the exchange of energy- 
momentum such that the total stress-energy tensor of the Universe 
is still conserved while the individual fluids feature an additional 
source term in their continuity equations: 



" M = ■ 



n)l(CDM) 



CT; 



(i(CDM) 



(1) 

(2) 



where C is a constant and <f> is the homogeneous value of the DE 
scalar field. These scenarios have been studied thoroughly concern- 
ing their effects on the background expansion history of the uni- 
verse ( |Amendola|2000l[2004l[Caldera-Cabral et al.|2009||Clemson| 
|et al.|[20TT|), their impact on the g rowth of linear density pertur- 



batio ns (|Arnendola & Quercellini||2003 Pettorino & Baccigalupi 
[20081 [Koyama et al.|2009||La Vacca et al.|2009| |Amendola et al 



|201 1| >, and their effects on the nonlinear dynamics of collapsed 
structures JMaccio et al.|2004||Bludi et al.|2010[[Li & Barrow|201 1[ 
|Baldi|2 01 la 2012c). Such investigations have provided relatively 
tight constraints on the maximum allowed values of the interaction 
constant C based on a wide range of observational probes (see e.g. 
|Bean et al.|2008l|Baldi & Pettorino|20TT]|Baldi & Viel|2010||Xra| 
[20091 ). 

It goes beyond the scope of the present paper to provide a 
comprehensive review of the basic properties and the main observa- 
tional constraints of coupled DE models, and we refer the interested 
reader to the mentioned literature for a more detailed discussion of 
these issues. Here instead we want to focus on the fact that present 
constraints seem to exclude with high statistical confidence the pos- 
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sibility of couplings with strength comparable to the gravitational 
interaction. If we define, following the literature (e.g. |Amendola| 
2000) the dimensionless coupling f) through the relation: 



C 



2 1 

3 Mp 



(3) 



where Mpi = 1/ V8tvG is the reduced Planck mass with G the 
Newton's constant, present observational bounds restrict the cou- 
pling to P < 0.15 at most (see e.g. |Bean et al.|2008[ |Xia|2009| 



|Baldi & Viel|2010||La Vacca et al.|2009) , corresponding to a scalar 



interaction constant 4G/3 /3 ~ 0.03G (see|Amendola 



rrespc 
|Amei 



2004) with 



a few percent of the gravitational strength 

All these observational constraints apply to standard coupled 
DE models, where the coupling C is a constant and where only 
one CDM particle species is present in the universe. More com- 
plex scenarios, featuring e.g. a time evolution of the coupling ( 
a "bouncing" DE self-interaction potential (as e.g. in 



Baldi 



Baldi 



2012a) or the existence of more than one CDM species interact- 



ing with DE iBrookfi eldet al.|20 08, Baldi 2012b) have been pro- 
posed in recent years and have been shown to easily evade such 
constraints, allowing for larger values of the coupling. In particu- 
lar, Brookfield et al. ( 2008 ) investigated a general framework where 
multiple CDM fluids interact with individual couplings to the DE 
scalar field, highlighting how the presence of multiple dark matter 
species determines new critical points of the background dynami- 
cal system for which the total effective background coupling is sup- 
pressed, thereby reducing the impact of the coupling on the cosmic 
expansion history. More recently, Baldi ( 20 1 2b ) studied in large de- 
tail a specific realization of such general framework, where CDM 
is composed by only two different particle species completely de- 
generate in all their physical properties except for the sign of their 
interaction constant with the DE scalar field. Such specific Multiple 
Dark Matter coupled DE model (or Multi-coupled DE, McDE here- 
after) is described at the background level by the following system 
of dynamic equations: 



■ + 3Hcj> - 



dV 



p+ + 3Hp+ 
p- + 3Hp- 

pr + 4:Hp r 

3H 2 



+Cp+ - Cp- 
-C<j>p+ , 

o, 

1 

M 



2 (Pr + P++ P-+ Pep) , 



(4) 

(5) 
(6) 
(7) 

(8) 



where an overdot represents a derivative with respect to the cosmic 
time t, H = a /a is the Hubble function, the two CDM species 
are denoted by the subscripts "+" and "— " according to the sign 
of their coupling constant, and the total CDM density is given by 

PCDM = P+ + P-- 

In the present paper, we will focus on this particular class of 
models which has the appealing feature of requiring no additional 
free parameters with respect to a standard coupled DE scenario 
with constant coupling. Without loss of generality for the purposes 
of our discussion, we will not include in the analysis the uncoupled 
baryonic component of the universe, and we will always assume an 
exponential form for the scalar self-interaction potential V(<f>) ap- 
pearing in the Klein-Gordon equation |4](. The background dynam- 
ics of this scenario has been described in detail in Baldi (2012b) 
which showed how for reasonable initial conditions (i.e. without a 
too large asymmetry between the two different CDM species in the 
early universe) the expansion history of the model is completely in- 
distinguishable from ACDM even for dimensionless couplings as 
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Figure 1. The evolution of the DE scalar field <j> as a function of the 
e-folding time In a in an uncoupled scenario (solid) and in McDE cos- 
mologies with different values of the coupling (3 (dotted, dashed, and dot- 
dashed). As one can see from the figure, after the initial evolution the scalar 
field gets trapped in the minimum of its effective self-interaction potential 
during most of matter domination (represented by the grey-shaded area) and 
during this period behaves like an uncoupled field. 



large as (3 — 10. This is due to the fact that during matter dom- 
ination the scalar field is trapped in the minimum of its effective 
potential V e e defined by: 

rfVcff = dV 
dcj> dcj> 



Cp+ + Cp- 



(9) 



and the scalar field velocity (f> is consequently very small, thereby 
suppressing the effect of the interaction in Eqs. \5\(s\ . 

In Fig. [Tj we show the evolution of the scalar field <f> as a 
function of the e-folding time In a for a range of McDE models 
with different couplings starting from the same initial conditions 
at high redshift. As the plot shows, even for coupling values 
larger than unity (dotted, dashed, and dot-dashed curves) the 
scalar field is attracted to a fixed point where it remains constant 
during the whole matter dominated epoch, identified in the plot 
by the grey-shaded area, thereby behaving like an uncoupled 
system (solid curve). In other words, the background evolution 
of the universe can be completely screened from arbitrarily large 
coupling values if two CDM species with an opposite sign of the 
coupling both interact with the DE scalar field. Such screening 
mechanism therefore does not allow to distinguish a McDE model 
from an uncoupled cosmology - or even from ACDM - through 
geometrical probes of the cosmic evolution. 



However, Brookfield et al. ( 2008 ) already showed that this de- 
generacy can be easily broken by the evolution of linear density 
perturbations. In fact, even for a vanishing effective coupling in 
the background equations Pj8j ), linear density perturbations fea- 
ture a significantly enhanced growth with respect to ACDM in 
the presence of such large coupling values. This is due to the fact 
that the symmetry between the two CDM species that holds in the 
background as a consequence of the new matter-dominated critical 
point, is broken at the level of density perturbations by any oscil- 
lation of the field around the minimum of the effective potential, 
which determines a different friction term for the fluctuations in 
the two different CDM species (see Baldi 2012b). This can be un- 
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derstood by having a look at the linear perturbations equations for 
the specific case of two CDM species with opposite couplings: 



-2H 



5- = -2H 



HV6 



S+ +4ttG [p-6-T R + p+S+T, 



(10) 



1 + 0- 



5- +^g\ P -S-Ta + P+ s+r R ] . (ii) 



In Eqs. {10|11| > the Y factors are defined as: 

r A = i + ^p 2 , r R = i-^p\ (12) 

and represent attractive (r^) or repulsive (F_r) corrections to grav- 
ity due to the long-range fifth-force mediated by the DE scalar field, 
while the second terms in the first squared brackets on the right- 
hand side are the friction terms associated to momentum conserva- 
tion, that break the symmetry of the two equations. From Eqs. ([To]- 
12 1 it is then clear that any dynamics of the scalar field (i.e. 0^0) 
will induce a different growth of the density fluctuations in the two 
CDM species, which for couplings larger than ft — V3/2 will 
then start repelling each other due to their mutual repulsive scalar 
force (see Baldi 2012b| for a more detailed discussion on the origin 
of the repulsion between fluctuations in the two different species). 
As a result, the perturbation in the more overdense species will 
keep growing while the other will start decaying, thereby produc- 
ing isocurvature modes even when starting from purely adiabatic 
initial conditions. In any case, the growth of overdensities in one 
of the CDM species will always be faster than the depletion of the 
corresponding underdensities in the other species generated by the 
repulsive scalar force, due to the additional attractive pull of stan- 
dard gravity. As a result, the total CDM density perturbations given 
by 



3CDM 



+ + H 



(13) 



will still grow in time with a rate that will significantly depend on 
the coupling strength and on the level of adiabaticity of the pertur- 
bations set. 

Such peculiar behavior is well described by Fig. [2] where we 
show the evolution of a gaussian density perturbation in real space 
in the context of four different McDE models with different cou- 
pling values. The plot displays the overdensity profile 8(x) as a 
function of the radial comoving coordinate x in arbitrary units. Dif- 
ferent colors correspond to different redshifts (green for z = 8, red 
for z = 2, black for z — 0) while dashed and solid lines cor- 
respond to the positively and negatively coupled species, respec- 
tively. The grey solid lines show the total CDM overdensity profile 
|T3j for the combined fluid, and the three redshifts are identified 
on these lines by different symbols (diamonds for z = 8, triangles 
for z — 2, and squares for z — 0). For all the different cases we 
start our integration at very high redshifts (z w I0 7 ) from an adia- 
batic configuration of the perturbations, such that £1+5+ — Q-S-. 
As the plot shows, in the absence of coupling (left panel) the over- 
density profiles in the two CDM species evolve in the same way, 
and the perturbations remain adiabatic until z — 0. Correspond- 
ingly, the total CDM overdensity follows the same evolution as the 
two individual species. On the contrary, in the presence of a non- 
vanishing coupling the evolutions of the two CDM species start to 
deviate from each other as time goes by. The negatively-coupled 
species starts to grow faster than its positively-coupled counter- 
part, due to the mutual repulsion given by the scalar fifth-force, 
thereby inducing the formation of isocurvature modes. In particu- 
lar, for the strongest coupling j3 — 3/2, we can clearly see in the 



rightmost panel of Fig. [2] that already at high redshifts (z — 8, 
green curves) the positively-coupled species shows no overdensity 
at all while the negatively-coupled one has a significantly more pro- 
nounced overdensity profile than in the uncoupled case. The sit- 
uation evolves further at later times (red and black curves) with 
the positively-coupled species quickly developing an underdensity 
(i.e. a "void") in correspondence to the large overdensity of the 
negatively-coupled one. 

It is however very interesting to notice that, although the en- 
hancement in the growth (or in the decay) rate of the individual 
CDM species dramatically increases with the coupling /3, the evo- 
lution of the total CDM density perturbation results only mildly af- 
fected by the fifth-force, as one can see by comparing the amplitude 
of the overdensity profiles at z = in the different models, indi- 
cated by the horizontal dotted lines. Unless for the most strongly 
coupled case of /3 = 3/2, which exhibits a slight increase of the 
total amplitude at z = 0, all the other models do not show any 
appreciable change in the evolution of the total CDM perturbation. 
This interesting result indicates that any probe of linear structure 
formation based on the evolution of the total gravitational poten- 
tial at linear scales will not be able to distinguish a McDE scenario 
with a coupling as large as /3 = y3/2 (i.e. a model with a scalar 
fifth-force twice as strong as gravity) from a standard ACDM cos- 
mology. In other words, despite the strikingly different evolution 
of the density fluctuations in the two CDM species and the appear- 
ance of isocurvature modes at late times, the superposition of the 
two CDM fluids would effectively hide the DE-CDM interaction 
also at the level of linear perturbations, even for couplings of order 
unity and larger, thereby allowing for scalar forces of order gravity 
in the dark sector without any directly observational footprint. As 
already discussed in|Baldi 



(2012b), this result shows that the claim 



by |Brookfield et aT7| 1 2008 1 that linear density perturbations allow 
to break the background degeneracy between McDE and ACDM 
models is actually valid only for relatively large coupling values. 

It is therefore natural to try to extend the investigation of 
McDE models to the nonlinear regime of structure formation, in 
order to study how their associated strong long-range scalar force 
between CDM particle pairs, be it either attractive or repulsive, 
might affect the dynamics of collapsed objects at small scales, and 
whether the nonlinear evolution of structures can provide some di- 
rect evidence of the presence of such dark interactions that would 
remain otherwise completely hidden both at the background and 
linear perturbations levels. To this end, we have run the first N-body 
simulations ever performed for this particular class of DE models, 
and we will discuss their main features and results in the next Sec- 
tions. 



3 THE SIMULATIONS 

In order to run our simulations, we have made use of the mod- 
ified version by Bal di et al.| {20 10 1 of the widely used parallel 
Tree-PM N-body code GADGET {Springel||2"00"5) that allows to 
self-consistently simulate the evolution of structure formation pro- 
cesses in the context of generalized coupled DE cosmologies by 
including all the relevant features of such models, from their spe- 
cific background evolution to the mass variation of CDM particles, 
to the effects of the fifth-force and of the extra friction acting on 
individual particles. Such modified code has been widely used in 
the past to investigate coupled DE cosmologies with constant cou- 
plings {Baldi et al.|2"0"T0l|Baldi & Viel|2010||Baldi|201 la) variable 
couplings (Baldi 201 lb||Baldi & Pettorino|201 1| and with bounc- 
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Figure 2. The evolution of a gaussian density profile in real space for the positively- and negatively-coupled CDM species (dashed and solid lines, respectively) 
and for the total CDM density perturbation (grey solid lines) at redshifts z = , 2 , 8 (black, red, and green colors, respectively). The plot shows the overdensity 
profile S(x) as a function of a radial comoving coordinate x in arbitrary units. The grey squares, triangles, and diamonds are overplotted to highlight the shape 
of the total CDM density perturbation <5q D m at the three different redshifts considered in the plot. The horizontal dotted grey lines show the total amplitude 
of the CDM density perturbations at z = 0. In the left panel, corresponding to an uncoupled case = 0, no difference can be seen in the time evolution of the 
density profiles of the various components, while when progressively increasing the coupling value, moving to the right panels, one can start noticing a slight 
difference in the evolution of the positively- and negatively-coupled CDM species at low redshifts already for = 1 . The effect becomes more evident (and 
appears at earlier times) for larger values of the coupling, and for the most extreme case of = 3/2 the two CDM species develop opposite density profiles 
by the present time, thereby showing the formation of isocurvature modes. Nevertheless, the total CDM perturbation amplitude is only mildly affected by this 
strikingly different evolution of the two individual CDM components. 



ing DE potentials (Baldi 2012a), and has b een recently e mployed 
for the development of the CoDEC^Project I Baldi 2012c 1 that pro- 
vides publicly available numerical data from simulations of a wide 
range of coupled DE cosmologies. Here we apply for the first time 
the same code to the case of McDE, which does not require any 
particular modification with respect to the original implementation 
since the code was developed in the first place to allow for multi- 
ple families of massive particles with individual couplings. There- 
fore, differently from previous applications where the two matter 
families were identified with baryons and CDM (with the coupling 
of the former always set to zero) here we discard the baryons by 
switching off any hydrodynamical treatment of the corresponding 
particle type, and we identify the two matter families with the two 
different CDM species, simply setting the respective couplings to 

With this approach, we have run six cosmological N-body 
simulations of intermediate resolution for six different values of 
the coupling (fi = , 1/2 , a/3/2 , 1 , ^/3/2 , 3 /2) with the aim to 
highlight the main qualitative features of McDE models in the non- 
linear regime. The simulations follow the evolution of 2 x 256 3 
CDM particles in a periodic box of 100 Mpc/h aside, with a mass 
resolution at the starting redshift Zi = 99 of m± — 2.24 x 10 9 
Mq /h. The gravitational softening has been set to e g = 10 kpc / h, 
corresponding to about 1/30-th of the mean interparticle separa- 
tion. Higher resolution simulations would have requested a consid- 
erably higher computational cost that is not justified by the needs of 
the present preliminary investigation of McDE cosmologies, which 
is mainly aimed at sampling the parameter space of the modes and 
highlighting its most prominent features, and will be performed in 
future works. 



www.marcobaldi.it/research/CoDECS 



Initial conditions have been set by rescaling the amplitude 
of the linear power spectrum obtained with the Boltzmann code 
CAMB < |Lewis et aTl|2000) for a ACDM universe with WMAP7 
cosmological parameters iKomatsu et al. 2011) between last scat- 
tering (z ~ 1100) and the starting redshift of the simulations with 
the specific growth factor computed for each model by numerically 
integrating Eqs. «|T0]TTJ. However, since the growth rate of McDE 
models is practically indistinguishable from that of ACDM at high 
redshifts, the initial conditions of all the six simulations are vir- 
tually identical. The random-phase realization of the linear power 
spectrum in the initial conditions adopts the same random seed for 
all the models, such that the main shape of the resulting large-scale 
structures will be identical in all the runs. 



4 RESULTS 

We now move to discuss the main outcomes of our numerical simu- 
lations, focusing on the qualitative behavior of McDE cosmologies 
in the nonlinear regime of structure formation. As we will see in the 
remaining of the paper, our results show that only in this regime it is 
possible to identify some specific footprints that could allow to de- 
tect and possibly constrain a long-range fifth-force of gravitational 
strength in the dark sector, which would otherwise be completely 
hidden - in these scenarios - both at the level of the background 
and of the linear perturbations evolution. 



4.1 Large-scale structures in Multi-coupled Dark Energy 
cosmologies 

We start the analysis of our simulations by investigating the quali- 
tative shape and evolution of large-scale structures forming from 
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Figure 3. The distribution of the positively- and negatively-coupled CDM particles (red and black points, respectively) in a slice of 100 X 100 X 2 Mpc//i 
for different values of the coupling j3. A higher-resolution color version of these plots is downloadable at: www.marcobaldi.it/data/MDM_figures.tar 
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almost identical initial conditions in the context of the differ- 
ent McDE models under investigation. In Figure [3] we plot the 
projected position of particles of the positively- and negatively- 
coupled CDM species (as red and black points, respectively) in a 
slice of the cosmological box of the different simulations at z — 0. 
The thickness of the slice is 2 Mpc/h and the positively-coupled 
particles are plotted first such that they might be covered by the 
subsequent plot of the negatively coupled particles. The different 
panels refer to different values of the coupling as explained in the 
legend of the figure. 

The left panel in the upper row displays the particle distribu- 
tion in the absence of coupling (J3 = 0) and shows the familiar 
filamentary structure of the cosmic web developing in any hier- 
archical structure formation scenario. In particular, the red points 
representing positively-coupled CDM particles are barely visible 
(except in voids) as their spatial distribution closely follows that of 
the negatively-coupled particles represented by the black points: as 
expected, since both particle types feel the same total force given 
only by the standard gravitational interaction, their distributions 
are indistinguishable from each other. Interestingly, the situation 
does not seem to change appreciably also for a coupling value of 
/3 = 1/2 (right panel in the upper row), corresponding to a scalar 
fifth-force with a strength equal to one third of standard gravity, for 
which the qualitative shape of the large-scale structures appears in- 
distinguishable from the uncoupled case. This comparison already 
represents a very significant result in itself: a scalar fifth-force with 
one third the strength of standard gravity does not seem to have 
any striking effect on the formation of large-scale cosmic structures 
even in the nonlinear regime probed with full N-body simulations. 
More detailed investigations focused on the internal structure of 
collapsed objects might then be necessary to highlight possible ob- 
servational footprints of scalar forces of this strength in the context 
of McDE scenarios. We remind again here that a scalar fifth-force 
with strength G/3 between CDM particles is ruled out at more than 
5cr (Bean et al.|2008l |Xia|2009l |Baldi & Viel|2010) for standard 
coupled DE scenarios with only one CDM particle type. 

If we now move to progressively larger couplings, however, 
some distinctive effects clearly start to appear, and the two CDM 
particle species do no longer look distributed in space in the same 
way. In the middle row of Fig. [5] we display the cases of/3 = V3/2 
(left) and p = 1 (right). For the former model - corresponding to a 
fifth-force with the same strength as standard gravity - the effect is 
still very mild and the overall shape of the cosmic web appears still 
identical to the uncoupled case, but a larger number of red points 
representing positively coupled particles starts to appear, showing 
that the two CDM families are no longer identically distributed over 
the whole simulation box. The qualitative effect of the DE-CDM 
interaction becomes much more evident for a coupling of f3 — 1, 
for which a clear shift between the two CDM species appears when 
looking at the edges of large and elongated structures (like e.g. cos- 
mic filaments) showing an evident ongoing process of separation 
between the two CDM particle types. 

Further increasing the coupling value, the bottom row of Fig. [3] 
shows in the left panel the particle distribution for f} — ^/3/2, 
where the development of two almost equally-shaped cosmic net- 
works shifted from each other by a distance of a few Megaparsecs 
becomes apparent. This development of "mirror" structures in the 
cosmic web represents the first prominent feature that characterizes 
McDE models and that potentially allows to distinguish them from 
a standard minimally coupled cosmology. Furthermore, it is possi- 
ble to notice already by eye that the distributions of the two individ- 
ual CDM species also start showing some deviation from the over- 



all shape of the CDM large-scale structure that was instead roughly 
preserved in the previous plots. In particular, for (3 — y3/2 - 
corresponding to a fifth-force twice as strong as gravity - one can 
clearly observe the fragmentation of large structures (like e.g. mas- 
sive halos or filaments) into smaller and more peaked independent 
objects, and a significant depletion of cosmic voids. Finally, for a 
fifth-force three times stronger than gravity, represented by the right 
plot of the bottom row of Fig. [3] (0 = 3/2), the two CDM particle 
species are distributed in a way that only vaguely retains the mem- 
ory of the initial shape of the cosmic web, with isolated and very 
concentrated objects composed almost uniquely by one of the two 
different types of CDM particles. Very interestingly, the original 
shape of the large-scale structure can still be traced by following 
the extended voids that form in the regions between two "mirror" 
structures made of different CDM species that have originated from 
the same "seed" structure and that have subsequently departed from 
each other as a consequence of the strong repulsive force between 
the two particle families. 

The visual inspection of the six panels of Fig. [3] then clearly 
shows how the nonlinear regime of structure formation can reveal 
distinctive features of McDE models that would be otherwise com- 
pletely indistinguishable from ACDM in their background and lin- 
ear perturbations evolution, and therefore represents a unique han- 
dle to detect and possibly constrain the long-range scalar forces that 
arise in the context of these scenarios. In particular, this first qual- 
itative analysis shows that McDE models determine, for coupling 
values larger than unity, the formation and the progressive separa- 
tion of "mirror" cosmic structures mainly composed by particles 
of one single CDM species, and the corresponding onset of halo 
fragmentation processes at late times, when large "seed" structures 
composed by a mixture of the two CDM families split into smaller 
objects mainly composed by one single particle type. However, our 
results also show that for couplings smaller than unity, which in- 
clude also the case of a scalar force with the same strength of stan- 
dard gravity, the impact of the fifth-force on the overall shape of 
large-scale structures is still very mild even at the nonlinear level. 

4.2 The nonlinear matter power spectrum 

We now move to a more quantitative analysis of the structures 
forming in the context of the different McDE models by computing 
the nonlinear matter power spectrum of the different CDM species 
and of the total CDM density field directly from our simulations. 
The power spectrum of these three different components is com- 
puted in Fourier space by determining their individual density fields 
on a cartesian grid with the same size of the PM grid used for the 
long-range force computation of the simulations (i.e. with 256 grid 
nodes on each side) through a Cloud-in-Cell mass assignment al- 
gorithm. The nonlinear power spectra computed according to such 
procedure allow us to investigate the effects of the DE-CDM in- 
teraction on the mass distribution up to the Nyquist frequency of 
the PM grid, corresponding to a wave number /cNy ~ 8/i/Mpc, 
i.e. down to length scales relevant for present and upcoming weak 
lensing surveys. 

In Figure [4] we show for each model under investigation the 
nonlinear matter power spectra of the positively- and negatively- 
coupled CDM species (P+(fc) and P_(fc), respectively plotted as 
red and black solid curves) as well as the power spectrum of the 
total CDM density field (PcDM(fc), blue solid curve). As a ref- 
erence, we also overplot on each panel the nonlinear power spec- 
trum P\(fc) of the standard ACDM cosmology as a green dashed 
line. The bottom plot of each panel shows in red the ratio be- 



© 201 1 RAS, MNRAS 000.[7]fiT| 



8 M. Baldi 



10000F 



1000 



100 



10 



1.5 

< 
Q. 

1 1.0 



0.5 



10000 



1000 



100 



P= 1/2 




10000T 



1000 : 



P = 3 1/2 /2 



100 : 



10 r 



z=0.00 




z=0.00 



pyp. 



1.5 



pyp. 



a 1.0 ~- 



- Pp-Hm/Pa 



Q_ - P /p 

„ r - CDM A 

O.o — 



0.1 



1.0 

k [ h Mpc 1 ] 



10.0 



0.1 



1.0 

k [ h Mpc 1 ] 



10.0 



10 



1.5 

< 

3 1.0 



0.5 



P=1 




10000 F 



1000 r 



100 - 



P = (3/2) 1 



z=0.00 




0.1 



1.0 

k [ h Mpc 1 ] 
10000 



1000 



10.0 




1.0 
k [ h Mpc 1 ] 



10.0 



100 



10 



1.5 



I 10 



0.5 



P = 3/2 




z=0.00 



pyp 




0.1 



1.0 

k [ h Mpc 1 ] 



10.0 



Figure 4. The upper plot of each panel displays the matter power spectrum of the positively- and negatively-coupled CDM species P+(k) and P_ (fc) and 
of the total CDM component PcdmW. plotted as red, black, and blue solid lines, respectively, and the matter power spectrum of the standard ACDM 
cosmology overplotted for visual reference as a green dashed line, all at z = 0. The bottom plots show in red the ratio of the power spectra of the two different 
CDM components P_ (k)/ P+ (A;) and in blue the ratio of the total CDM power spectrum to the ACDM reference case, Pcdm / P\ ■ As the different panels 
show, for a coupling of j3 = 1 / 2 no effect can be seen in the plots, while for larger coupling values a progressive suppression of power at small scales start to 
appear, which develops an almost-constant plateau at scales k > 1 h/Mpc for ft ^ 1. Such plateau progressively moves back to the standard ACDM case for 
increasing couplings, and leaves a very characteristic maximum of the power suppression at intermediate scales between k ~ lh /Mpc and k ~ 4 h /Mpc. 
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tween the power spectra of the two different CDM species, i.e. 
P-(k)/ P+(k), while in blue it displays the ratio between the to- 
tal CDM power spectrum of each specific model and the reference 
ACDM case, P C DM(fc)/-P A (&). 

As one can see from the plots, for a coupling /3 = 1/2 the 
power spectra of the two different CDM species are practically in- 
distinguishable from each other and from the total CDM one, and 
no difference can be seen also when comparing the latter to the 
standard ACDM power spectrum. This result reinforces the conclu- 
sion that, as suggested by the visual inspection of the large-scale 
structures discussed above, a scalar fifth-force with one third the 
strength of standard gravity does not affect in any appreciable way 
the density field in the context of McDE cosmologies down to the 
scales that we are able to resolve with present simulations. 

For larger coupling values, however, some effects start to 
clearly appear both in the relative evolution of the matter power 
spectra of the two different CDM species, and in the the total CDM 
power spectrum. In particular, for P = V3/2, corresponding to a 
scalar fifth-force with gravitational strength, a suppression of the 
total CDM power at scales smaller than k ~ 1/i/Mpc clearly 
appears in the plot, and progressively increases for smaller and 
smaller scales up to a maximum suppression of ~ 40% at the 
Nyquist frequency, while for k > 1 /i/Mpc the effect is basically 
absent and all the different power spectra are still identical. This 
is a particularly interesting result since it highlights a direct ob- 
servational feature of a long-range scalar force with gravitational 
strength that appears only at scales that will be efficiently probed 
by the next generation of weak lensing surveys, while leaving larger 
scales as well as the background evolution of the universe com- 
pletely unaffected. Another interesting feature of this plot is the 
different behavior at small scales of the two CDM species, that 
confirms the evolution displayed in Fig. [2] based on a purely lin- 
ear treatment. More specifically, the plot shows how, starting from 
the smallest scales, perturbations that are initially adiabatic tend to 
develop some isocurvature component as a consequence of their 
different gravitational interactions determined by the friction terms 
in Eqs. pOfTT) . 

When moving to larger coupling values these effects become 
much more prominent, with a suppression of the total CDM power 
that can reach a level of about 40 — 50% already at k ~ 2 /i/Mpc 
while still showing no suppression at all at the largest scales tested 
with our simulations. It is nevertheless particularly interesting to 
notice that for couplings /3 ^ 1, the small-scale suppression of the 
matter power spectrum features an almost constant plateau whose 
level tends to move back to the standard ACDM amplitude for in- 
creasing values of the coupling, leaving a clear maximum in the 
power spectrum suppression at intermediate scales. Such maximum 
(i.e. a minimum in the ratio Pcdm /Pa in Fig.[4j also moves to pro- 
gressively larger scales for increasing values of the coupling such 
that for the most extreme case addressed in this work, ft = 3/2, 
the power spectrum is suppressed at the smallest scales by only 
about 10% and it still converges to the ACDM amplitude at the 
largest scales, while a maximum suppression of about 60% ap- 
pears at k ~ 0.4/i/Mpc. This result is particularly interesting 
since it identifies another very distinctive observational footprint 
of McDE models, and shows an effect that could have only been 
tested through full N-body simulations: the fragmentation of large 
collapsed objects made by a mixture of the two CDM particle 
types into smaller structures composed primarily by a single par- 
ticle type. 

This fragmentation is triggered by local deviations from per- 
fect spherical symmetry of the density profiles of the two different 



CDM species within individual initial "seed" overdensities. Such 
asymmetries arise as a consequence of tidal effects and merging 
processes, and then cannot be properly taken into account in the 
context of linear perturbation theory or even with a nonlinear spher- 
ical collapse treatment. In fact, if the density profiles of the two 
CDM species within a common "seed" overdensity are not both 
perfectly spherically symmetric, the two species will start repelling 
each other away also along some radial direction, thereby triggering 
the fragmentation of the "seed" structure into two smaller struc- 
tures. Such fragmentation process then determines an initial sup- 
pression of the nonlinear matter power, starting from the smallest 
scales which are the first ones to collapse, since the development of 
isocurvature modes and the separation of the sibling structures orig- 
inated from an initial common "seed" reduce the total amplitude of 
density perturbations. However, when two sibling halos which have 
originated from the same "seed" structure are finally separated, they 
are no longer subject to the repulsive force given by the mixture of a 
comparable amount of particles of the two types, and only feel the 
strongly enhanced self-attraction given by the sum of their stan- 
dard self-gravity and their attractive fifth-force, thereby growing 
significantly faster than in the context of standard gravity and pro- 
gressively becoming more overdense. This makes the small-scale 
power to grow faster than in ACDM in the "post-fragmentation" 
phase, thereby progressively reducing the initial suppression. Such 
process of initial suppression and subsequent faster growth of the 
power spectrum amplitude is expected to occur at larger and larger 
scales as time goes by, which also explains why the maximum sup- 
pression scale moves towards lower wavenumbers for increasing 
couplings. 

All these different peculiar features at scales potentially 
testable through galaxy clustering and weak lensing observations, 
and identified here for the first time, represent the only characteri- 
stic footprints of McDE cosmological models which might allow 
to directly constrain this type of scenarios and to possibly rule out 
long-range scalar forces in the dark sector, at least for sufficiently 
large interaction strengths. 



5 CONCLUSIONS 

In this work we have presented the results of the first N-body 
simulations ever performed for Multi-coupled Dark Energy cos- 
mological models, characterized by the existence of two distinct 
CDM particle species interacting with opposite couplings with a 
classical Dark Energy scalar field. Such models strongly alleviate 
the impact that standard coupled DE scenarios have on both the 
background expansion history of the universe and the growth rate 
of linear density perturbations, and result completely indistinguish- 
able from the standard ACDM cosmology over a large range of 
coupling values. As a consequence, long-range scalar fifth-forces 
of gravitational strength that would be starkly incompatible with 
basic cosmological and astrophysical observations within standard 
coupled DE scenarios, are instead fully viable at the background 
and linear perturbations level in the context of Multi-coupled DE 
models. 

We have therefore investigated for the first time the effects of 
these cosmologies in the nonlinear regime of structure formation, 
by means of a series of intermediate-resolution N-body simulations 
that implement all the relevant effects characterizing Multi-coupled 
DE models, with the aim to qualitatively highlight possible specific 
footprints of such scenarios and of their related long-range scalar 
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forces over a wide range of parameter values. In particular, we have 
studied the qualitative evolution of the large-scale structure shape 
and of the relative spatial distribution of the two different CDM 
families for different values of the DE-CDM coupling constant, 
and more quantitatively determined the impact of the different 
models on the nonlinear matter power spectrum at scales that will 
be efficiently probed by present and future weak lensing surveys. 

With the former analysis, we have shown how a scalar 
fifth-force with gravitational strength (or weaker) does not ap- 
preciably affect the overall shape of the large-scale structures 
forming in our simulated cosmological volumes, thereby making 
an extremely difficult task to distinguish Multi-coupled DE models 
from minimally coupled scalar field cosmologies and even from 
the standard ACDM model. Higher-resolution simulations will 
be necessary in order to study the impact of such models on the 
inner regions of massive collapsed halos. With larger values of 
the coupling, instead, some peculiar features of the interaction 
can be identified in the nonlinear regime already at the present 
resolution level, for scenarios that would still appear completely 
indistinguishable from ACDM in the background and in the growth 
of linear perturbations. In particular, our simulations show that in 
Multi-coupeld DE models the spatial distributions of the two dif- 
ferent CDM species start to deviate from each other for sufficiently 
large coupling values, thereby leading to the progressive separation 
of the two CDM families into two "mirror" cosmic webs of 
collapsed halos and filaments that progressively move away from 
each other. This peculiar process also leads to the fragmentation 
of large collapsed objects made of an almost equal mixture of the 
two CDM particle types into smaller structures mainly composed 
by one particle species only. Such effects become dramatic for 
fifth-forces more than three times stronger than standard gravity, 
and destroy the overall shape of the cosmic large-scale structures, 
thereby providing a direct way to constrain scenarios that would 
be otherwise indistinguishable from ACDM. 

For the latter and more quantitative analysis, we have com- 
puted the nonlinear matter power spectra of the two different CDM 
families and of the total CDM fluid directly from our simulations, 
and shown that a potentially measurable suppression of power at 
scales k > 1 h/Mpc starts to appear already for scalar fifth-forces 
with the same strength as standard gravity, and progressively 
becomes more significant for larger values of the DE-CDM 
coupling. Finally, although such small-scale suppression tends to 
decrease again for even larger coupling values, as a consequence of 
the enhanced growth of the separate "mirror" halos once they have 
departed from each other, a characteristic maximum suppression at 
intermediate scales clearly remains imprinted in the matter power 
spectrum and represents a very characteristic footprint that could 
allow to rule out such extremely large coupling values using future 
weak lensing surveys. 

We can therefore conclude that McDE models - which are 
completely indistinguishable from the standard ACDM cosmology 
both at the background and at the linear perturbations level for a 
wide range of coupling values - can have interesting and pecu- 
liar footprints in the nonlinear regime of structure formation. How- 
ever, such effects become particularly prominent, and might pos- 
sibly rule out or constrain the model, only for relatively large val- 
ues of the dimensionless coupling, while for scalar interactions up 
to gravitational strength they remain sufficiently subtle to evade 
all the possible diagnostics that are allowed by the resolution of 



present simulations. A more detailed investigation of these models 
is therefore required to explore their possible effects on the strongly 
nonlinear systems represented by the internal regions of massive 
collapsed halos, and to provide a more quantitative determination 
of their observational signatures. 
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